## This program generates Figure 6 from "Financial Constraints, Sectoral Heterogeneity,
##  and the Cyclicality of Investment" by Cooper Howes

##### Step 0: Preliminaries #####

## Clear variables and data
rm(list=ls())

## Load packages
library(tidyverse)

##### Step 1: Load data #####

## Load impulse responses from data (created in program to generate Figure 2)
qfr_irfs <- read_csv("irfs/qfr_irfs.csv")

## Load impulse responses from model (created in Dynare file "baseline_model.mod")
irfdata <- t(read.csv("irfs/shock_irfs_baseline.csv",header=F))

## Extract names and relabel
colnames(irfdata) <- irfdata[1,]
rownames(irfdata) <- NULL
irfdata <- irfdata[-1,]

##### Step 2: Plot figures #####

## Create x-axis labels
qtrs <- 0:16

## Specify y-axis range
yrange <- c(-0.02,0.045)

## Create new window to hold both panels
windows(width=12,height=5)
par(oma=c(5,1,0,1))
m <- matrix(1:2,nrow = 1)
layout(mat = m)

## Plot IRFs for capital stocks in model
par(mar = c(2,2.5,3,1))
plot(x=qtrs,y=irfdata[,"k_d_e_m"],col="blue",type="b",pch=1,lwd=4,ylim=yrange,ann=FALSE,axes=FALSE,frame=TRUE)
lines(x=qtrs,y=irfdata[,"k_n_e_m"],type="b",pch=2,col="red",lwd=4)
lines(x=qtrs,y=irfdata[,"k_e_m"],type="b",pch=3,col="orange",lwd=4)
abline(a=0,b=0,col="black",lwd=1)
axis(1,cex.axis=1.25)
axis(2,cex.axis=1.25,at=pretty(yrange),lab=paste0(pretty(yrange)*100, "%"),las=TRUE)
title(main="Model capital stock", col.main="blue", font.main=2,cex.main=1.5,line=1)
mtext("Quarters after shock", side=1, line=2.5,cex=1.25)

## Plot IRFs for capital stocks from the data
par(mar = c(2,2.5,3,1))
plot(x=qtrs,y=qfr_irfs$dur,col="blue",type="b",pch=1,lwd=4,ylim=yrange,ann=FALSE,axes=FALSE,frame=TRUE)
lines(x=qtrs,y=qfr_irfs$nondur,type="b",pch=2,col="red",lwd=4)
lines(x=qtrs,y=qfr_irfs$total,type="b",pch=3,col="orange",lwd=4)
abline(a=0,b=0,col="black",lwd=1)
axis(1,cex.axis=1.25)
axis(2,cex.axis=1.25,at=pretty(yrange),lab=paste0(pretty(yrange)*100, "%"),las=TRUE)
title(main="Data capital stock", col.main="blue", font.main=2,cex.main=1.5,line=1)
mtext("Quarters after shock", side=1, line=2.5,cex=1.25)

## Add legend at the bottom
par(fig = c(0, 1, 0, 1), oma = c(0, 0, 0, 0), mar = c(0, 0, 0, 0), new = TRUE)
plot(0, 0, type = "n", bty = "n", xaxt = "n", yaxt = "n")
legend(x="bottom",inset=0,legend=c("Durable","Nondurable","Total"),
       lwd=2,lty=1,pch=c(1,2,3),col=c("blue","red","orange"),cex=1.5,horiz=TRUE)


